Estimation of initial conditions from a scalar time series 
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We introduce a method to estimate the initial conditions of a mutivariable dynamical system from 
a scalar signal. The method is based on a modified multidimensional Newton-Raphson method which 
includes the time evolution of the system. The method can estimate initial conditions of periodic 
and chaotic systems and the required length of scalar signal is very small. Also, the method works 
even when the conditional Lyapunov exponent is positive. An important application of our method 
is that synchronization of two chaotic signals using a scalar signal becomes trivial and instantaneous. 
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A trajectory that a given dynamical system traverses 
in its state space depends on the particular set of ini- 
tial conditions with which it starts. In particular, the 
state of a chaotic system at a latter time is exponen- 
tially sensitive to changes in its initial state 0. This 
defining feature of a chaotic system leads to a complex 
behaviour in state space that appears random yet is de- 
terministic, which means that an initial state uniquely 
fixes the future course of its evolution. Though their are 
several invariant measures of a chaotic system which are 
not sensitive to the intial conditions, the exact trajectory 
crucially depends on the intial state and hence is difficult 
to reproduce due to sensitivity to initial conditions. 

In light of these facts, it is interesting and impor- 
tant to ask whether the complete set of initial conditions 
of a given multivariable dynamical system can be esti- 
mated from a given scalar time series for a single state 
space variable. We show that this question can be an- 
swered in the affirmative and present a novel and simple 
method to estimate the initial conditions. Our method is 
based on a modified multidimensional Newton-Raphson 
method |l|,|| where we include the time evolution of the 
system. The length of the time series required for the 
calculations is typically very small. 

Our results raise some interesting issues regarding the 
information content of a time series. In standard em- 
bedding techniques, a vector space is constructed from 
successive iterates of a single variable and a trajectory is 
reconstructed in this space ||. While embedding, it is 
crucial to choose an appropriate time delay so that the 
successive iterates are well resolved and contain quali- 
tatively different information. In our method we use a 
very small time series and the total duration is typically 
much less than the standard delay time in embedding 
techniques. It is interesting that we can recover the ini- 
tial conditions and hence the trajectory from such a small 



duration of time series. 

An important application of our method is in the prob- 
lem of synchronization of two chaotic systems Q . This 
problem itself has attracted wide attention in recent 
times due to its potential application to secured com- 
munication and parameter estimation |^||. Our 
estimation of the initial conditions makes the problem 
of synchronization almost trivial. We also find that our 
method works for most of the cases where other methods 
fail |0. 

Let us consider an autonomous dynamical system 
given by, 

x = F(x), (I) 

where x = (x\,X2, ■ ■ ■ ,x c i) is a d-dimensional state vec- 
tor whose evolution is governed by the function F = 
(Fi, i<2, . . . , Fd). Given an initial state vector x(0) at 
time t = 0, the time evolution x(i) is uniquely deter- 
mined by Eq. (Q) . Now let us assume that only one com- 
ponent of the state vector is known to us and we take it 
to be x\(t) without loss of generality. The problem that 
we address is to obtain the initial state vector x(0) from 
the knowledge of the scalar signal X\{t). 

Let y(0) denote a random initial state vector and y(t) 
its time evolution obtained from Eq. (Q). Let w(t) denote 
the difference 

w(t) = y(t) - x(t). (2) 

We look for the solution of the equation 

w(t) = 0. (3) 

Noting that the initial state vectors y(0) and x(0) 
uniquely determine the difference w(i), one of the solu- 
tions of Eq. (|^) is y(0) — x(0) = and this is the solution 
that we are searching for. 
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We now introduce the notation w n = w™(y°,x°) = 
w(nAt), where At is a small time interval. Similarly, 
y" = y{nAt) and x" = x(nAt). With this notation 
condition (||) can be written as w" = 0. 

Our approach to the solution of Eq. (||) is a modified 
Newton- Raphson method Q which includes the time evo- 
lution of the system. 

Let us first consider w 1 . We have 

= w 1 (x°,x°), 
= w 1 (y + <5y°,x ), 

= wV,x°) + (5y°.V y o)w 1 (y ,x ) + O((Sy ) 2 ), (4) 

where <5y° = x° — y° = — w° and the last step is a Taylor 
series expansion in 5y . For small At, we can write 

w^yV ) =w° + Ai[F(y°)-F(x )] + O((Ai) 2 ). (5) 

Substituting Eq. (B) in Eq. (Q) and neglecting higher or- 
der terms, we get, 



wV.x ) = w° - At(w°.V y o)F(y°) 



(6) 



It is convenient to write the above equation in a matrix 
form as 



W 1 = {I + AtJ°)W°, 
= A°W°, 



(7) 



where W n is the column matrix corresponding to the vec- 
tor w n , I is the identity matrix, A n — I + AtJ n , and 
the elements of the Jacobian matrix J n are given by 



_ dF t (y n ) 



(8) 



Next we consider w 2 or W 2 . Proceeding as above, we 
get (see Eq. (§), 

W 2 = {I + AtJ^W 1 

= (I + At J x )(/ + At J°)W° 

= A 1 A°W°. (9) 

Similarly, the equation for W n is 

W n = {I + AtJ n - 1 )(I + AtJ n - 2 )---(I + AtJ°)W° 
= A n ~ 1 A n ~ 2 ■ ■ ■ A°W°. (10) 

We now concentrate on the first component of the sig- 
nal whose time series is assumed to be known. For a 
c?-dimcnsional system we need d — 1 equations to deter- 
mine the initial state vector x°. Eqs. (0), (||) and jic| ) 
give us the required relations. 



W 2 = ^ u A%Wl 
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These are d — 1 simultaneous equations for W°. 

The numerical procedure is as follows. We set the in- 
tial state of system (|l|) to a random initial guess vector 

(y°)oid with (y°)<>id = x i and evolve it: usin s E q- d)- 

Using this vector y(t) we write down d — 1 simultaneous 
equations (Eqs. (|Tl|)) which can be solved for d — 1 un- 
known components of w° = — Sy°. Also, Sy® = 0. Thus 
the initial guess vector can be improved by 



(A 



' old 



Sy°. 



(12) 



This sets up an iterative scheme giving us better and 
better estimates of the initial vector which converge to 
x°. 

We note that as in Newton-Raphson method, the 
choice of the initial guess vector can be very impor- 
tant In some cases, the iterative procedure of Eq. ( p^ ) 
may not converge or converge to a wrong root. In such 
cases, a different choice of initial guess vector can be use- 
ful. 

We further note the similarity of our method with the 
so called method of variational equations in analytical 
dynamics |]ll|] . The method of variational equations can 
be appied to a known Hamiltonian system to determine 
an unknown neighbouring trajectory to an already known 
one. There, the method requires a complete particular 
solution of a known set of Hamiltonian equations of mo- 
tion. In contrast, we have used our method for dissipative 
chaotic systems. In such systems an analytical solution 
of the equations of motion cannot be known. Further our 
method requires only one component of a complete tra- 
jectory to be sampled. This has important consequences 
in the problem of synchronization using a scalar signal. 

We now demonstrate our method of estimating the ini- 
tial state. As our first example we discuss the Rossler 
system given by Q , 

x\ = -x-x ~ x 3 , 
±2 = X\ + ax 2 , 

x 3 =b + x 3 (xi-c). (13) 

First, we consider a case when the time series for x\ 
is given and we want to estimate (x^xjj). We chose the 
parameters (a, 6, c) such that the system is in the chaotic 
regime and the initial state x° is on the chaotic attractor. 
We start with an arbitrary initial state y° = (yj, 2/2,2/3) 
with 2/1 = x\. From Eqs. ( pi] ) we get a pair of simulta- 
neous equations as, 
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w\ = y\ -x\ = At6y% + AtSy° 
w \=yl-xl= (2At + a{At) 2 ) 5yl 

+ (2At+(y° 1 -c)(At) 2 )6y 3 (14) 

which can be solved for (Sy^dy®). With 8y® — we 
use these in an iterative manner (Eq. (|l2|)) to obtain the 
correct intial conditions. 

Tabic 1 shows successively corrected (y 2,1/3) obtained 
using the iterative process as discussed. These are the 
successive estimates for (x^x®). Let = \y® — xf\ denote 
the absolute error in the estimation of x®. In Fig. 1(a) we 
plot a graph corresponding to Table 1 showing errors 
and e-i (on logarithmic scale) plotted against the number 
of iterations of our method (Eq. (|l2|)). From Table 1 and 
Fig. 1(a) we see that the successive estimates converge to 
the correct values of (x®, x®). Using only two data points 
in the given time series x\(t), we can thus readily esti- 
mate the full initial state x°. We also note that the rate 
of convergence is very good. In about 8 to 10 iterates we 
obtain the initial values (x^x®) to within computer ac- 
curacy. If we write the deviations of the succesive iterates 
from the correct values in the form 

(e i ) n = \(y°) n -x i \~ e -™, (15) 

where n is the number of iterations, then the value of the 
parameter a is found to be 2.01 for e2 and 2.02 for e^. 
This is consistant with the fact that Newton-Raphson 
method has a quadratic convergence jl],^). 

We note that the largest Lyapunov exponent for the 
subsystem (y§,y§) (conditional or subsystem Lyapunov 
exponent) is positive [Q . The success of our method does 
not depend on whether this Lyapunov exponent is pos- 
itive or negative. This is important for synchronization 
of chaotic signals as will be discussed afterwards. 

We next present cases where time series for the vari- 
ables X2 and X3 of the Rossler system are given. The 
procedure is similar to the case of time series for x\ as 
discussed above. Fig. 1 (b) shows the errors e\ and e 3l 
when time series for X2 is given, plotted a gain st the num- 
ber of iterations. The parameter a (Eq. (fig)) is 1.97 for 
ei and 1.95 for e 3 . This again indicates a quadratic con- 
vergence. Similarly, Fig. 1 (c) shows the quantities e\ 
and e2 when time series for X3 is given, as a function of 
the number of iterations. The parameter a (Eq. (|T^)) 
is 1.28 for ei and 1.30 for e2 which shows a convergence 
slower than quadratic. We note that the largest sub- 
system Lyapunov exponent is negative when time series 
for X2 is given and is positive when time series for X3 is 
given [Q. 

As our next example we consider the Chua's circuit 
which in its dimensionless form is given by |l3| ] , 

±1 = a(x 2 - xi - f{xi)), 

x 2 = Xl - X 2 + X3, 
±3 = -f3x 2 , 



f( Xl ) = bx x + i(a - b) [\xi - 1| - \xi + 1|] . (16) 

We chose the parameters a, b, a and (3 such that the at- 
tractor is a limit cycle. The initial state x° is chosen in 
the basin of attraction of this limit cycle. Fig. 2 shows 
the errors e2 and e3 when a time series for the variable x\ 
is given as a function of the number of iterations of our 
method. The parameter a (Eq. (|l5|)) in this case takes 
values 1.29 for e2 and 1.25 for e3 showing a slower than 
quadratic convergence. 

We have also applied our method for the cases when 
time series for the variables X2 and X3 are given and also 
for cases when parameters are such that the attractor is 
chaotic. In all the cases we are able to estimate the full 
initial state vector. 

We have successfully applied our method to estimate 
the initial state vector using a given scalar time series 
for many other dynamical systems as well. These in- 
clude Lorenz system p"4| in its periodic, chaotic or in- 
termittent regimes, the disk dynamo system modelling a 
periodic reversal of earth's magnetic field |l5|,[l6|], a 3-d 
plasma system formed by a three wave resonant coupling 
equations |]l7f and a four dimensional phase converter 
circuit [pfl. 

Now we will discuss an important application of our 
estimation method in the problem of synchronization of 
two identical chaotic systems coupled unidirectionally by 
a scalar signal. Let us suppose that Eq. (Q) describes a 
chaotic system and let us consider a replica of it given by, 
y = F(y). Without losing generality let's further assume 
that a scalar output signal x\{t) is given. The aim is to 
synchronize vector y(t) with x(t) using this scalar signal. 
Our method of estimating initial conditions makes this 
procedure trivial. Using the scalar signal we estimate the 
initial vector x° and set y° = x°. This clearly leads to 
an instantaneous synchronization of the two trajectories. 
As we have demonstrated in the case of Rossler system, 
our method works even when the largest conditional Lya- 
punov exponent is positive which is the case when other 
methods of synchronization are known to fail [§,[l0| . 

To summarise, we have introduced a novel yet simple 
method to estimate initial conditions of a multivariable 
dynamical system from a given scalar signal. Our method 
is based on a multidimensional Newton-Raphson method 
where we include the time evolution of the system. The 
method gives a reasonably fast convergence to the cor- 
rect initial state. The required length of the time series 
is very small. The method works even when the largest 
conditional Lyapunov exponent is positive. An impor- 
tant consequence of the method is that the problem of 
synchronization of identical chaotic systems using scalar 
signal becomes trivial since evolution of two such systems 
can then be started from identical initial states. 
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TABLE I. The table presents successive estimates of the unknown components of the initial state vector, i.e. (j/JJ, 3/3), in a 
Rossler system for which a time series xi(t) is given, n is the number of iterations of our method (Eq. (^)). A clear convergence 
is seen towards the actual values x%, x%, written at the bottom of the table. The parameters are (a, b, c) = (0.2, 0.2, 9.0). The 
time step At is 0.01 and we use fourth order Runge-Kutta method for time evolution of Eq. (|l3|). The corresponding plots are 
shown in Fig. 1(a). 
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X2 = -5.178173232007801 


xl = 0.4759645128337372 
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Fig.1 (AM and REA) 

FIG. 1. Plot (a) shows the errors ei and e$ on a logarithmic scale, as a function of n, the number of iterations of our method 
for Rossler system when time series for variable xi is given. The crosses show the values for e2 and solid triangles those for 
ez. Errors are seen to approach zero as n increases. The values are taken from Table 1. Similarly plot (b) shows the errors ei 
(crosses), and ez (solid triangles) as a function of n when the time series for X2 is given. Plot (c) shows the errors e\ (crosses), 
and e2 (solid triangles) as a function of n when the time series for £3 is given. For all these cases, the parameters are such that 
the attractor is chaotic and the initial state is on the attractor. 
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Fig.2 (AM and REA) 



FIG. 2. The figure shows errors e2 and e3 on a logarithmic scale, as a function of n, the number of iterations of our method 
for Chua's circuit when time series for the variable x\ is given. The crosses show the values for e2 and solid triangles those 
for ez. Errors are seen to approach zero as n increases. The parameters are(4, 3, 2, 0.1758) and are such that the attractor is a 
limit cycle and the initial state is in the basin of the attractor. 
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